PyTorch 单车预测器

共享单车公开数据集（Capital Bikeshare）。

该数据是从2011年1月1日到2012年12月31日之间某地的单车使用情况，每一行都代表一条数据记录，共17 379条。一条数据记录了一个小时内某地的星期几、是否是假期、天气和风速等情况，以及该地区的单车使用量（用cnt变量表示），它是我们最关心的量。

Note

本实例来自《深度学习原理与PyTorch实战（第2版）》，本文是我对该节的自学实践。

预测接下来一段时间，该地区单车数量的走势。

导入相关包：

import numpy as np
import pandas as pd # 读取CSV文件的库
import torch
import torch.optim as optim
import matplotlib.pyplot as plt
# 直接在Notebook中显示输出图像
%matplotlib inline

科学家早已从理论上证明，用有限多的隐含神经元可以逼近任意的有限区间内的曲线，这叫作通用逼近定理（universal approximation theorem）。

数据集

数据集位于这里，下载 data.csv。

加载数据集：

# 读取数据到内存，rides为一个dataframe对象
data_path = 'drive/MyDrive/hour.csv' 
rides = pd.read_csv(data_path)
rides.head() # 输出部分数据
counts = rides['cnt'][:50] # 截取数据
x = np.arange(len(counts)) # 获取变量x
y = np.array(counts) # 单车数量为y
plt.figure(figsize = (10, 7)) # 设定绘图窗口大小
plt.plot(x, y, 'o-') # 绘制原始数据
plt.xlabel('X') # 更改坐标轴标注
plt.ylabel('Y') # 更改坐标轴标注

上面代码将会绘制出数据集中的前五十个数据，即不同时间点下的单车数量：

损失函数

对于所有的数据样本，神经网络预测的单车数量与实际数据中单车数量之差的平方和的均值：

L = \frac{1}{N} \sum_{i = 1}^{N} ({\hat{y}}_{i} - y_{i})^{2}

这是一个计算平均平方误差（Mean Squared Error, MSE）的公式，常用于回归问题中评估模型的性能。

单车预测器1.0

构建一个单一隐含单元的神经网络，采用 10 个隐含单元。网络搭建：

# 输入变量，1,2,3,...这样的一维数组
x = torch.FloatTensor(np.arange(len(counts), dtype = float))

# 输出变量，它是从数据counts中读取的每一时刻的单车数，共50个数据点的一维数组，作为标准答案
y = torch.FloatTensor(np.array(counts, dtype = float))

# 设置隐含神经元的数量
sz = 10 

# 初始化输入层到隐含层的权重矩阵，它的尺寸是(1,10)
weights = torch.randn((1, sz), requires_grad = True)

# 初始化隐含层节点的偏置向量，它是尺寸为10的一维向量
biases = torch.randn((sz), requires_grad = True)

# 初始化从隐含层到输出层的权重矩阵，它的尺寸是(10,1)
weights2 = torch.randn((sz, 1), requires_grad = True)

迭代时对网络进行训练：

# 设置学习率
learning_rate = 0.001 
# 该数组记录每一次迭代的损失函数值，以方便后续绘图
losses = []

x = x.view(50,-1)
y = y.view(50,-1)
for i in range(100000):
	# 从输入层到隐含层的计算
	hidden = x * weights + biases
	
	# 此时，hidden变量的尺寸是(50,10)，即50个数据点，10个隐含神经元
	# 将sigmoid函数应用在隐含层的每一个神经元上
	hidden = torch.sigmoid(hidden)
	
	# 隐含层输出到输出层，计算得到最终预测值
	predictions = hidden.mm(weights2)
	
	# 此时，predictions的尺寸为(50,1)，即50个数据点的预测值
	# 通过与数据中的标准答案y做比较，计算均方误差
	loss = torch.mean((predictions - y) ** 2)

	# 此时，loss为一个标量，即一个数
	losses.append(loss.data.numpy())
	
	if i % 10000 == 0: # 每隔10 000个周期打印一下损失函数数值
		print('loss:', loss)
	
	# *****************************************
	# 接下来开始执行梯度下降算法，将误差反向传播
	loss.backward() # 对损失函数进行梯度反传
	
	# 利用上一步计算中得到的weights、biases等梯度信息更新weights和biases的数值
	weights.data.add_(- learning_rate * weights.grad.data)
	biases.data.add_(- learning_rate * biases.grad.data)
	weights2.data.add_(- learning_rate * weights2.grad.data)
	
	# 清空所有变量的梯度值
	weights.grad.data.zero_()
	biases.grad.data.zero_()
	weights2.grad.data.zero_()

执行后，可以看到随着训练，损失减小的过程：

loss: tensor(2152.5972, grad_fn=<MeanBackward0>)
loss: tensor(668.6163, grad_fn=<MeanBackward0>)
loss: tensor(532.7995, grad_fn=<MeanBackward0>)
loss: tensor(492.9125, grad_fn=<MeanBackward0>)
loss: tensor(474.9532, grad_fn=<MeanBackward0>)
loss: tensor(471.4581, grad_fn=<MeanBackward0>)
loss: tensor(470.8375, grad_fn=<MeanBackward0>)
loss: tensor(470.4477, grad_fn=<MeanBackward0>)
loss: tensor(470.0748, grad_fn=<MeanBackward0>)
loss: tensor(469.7846, grad_fn=<MeanBackward0>)

其中，我们进行了十万步的训练迭代，在每一次训练中，都使用数据集的前 50 个数据。

绘制损失曲线：

plt.plot(losses)
plt.xlabel('Epoch')
plt.ylabel('Loss')

注：我的图与书上略有不同，两者趋势一致，但书上是一条细线，我这里感觉就跟钢笔漏水了一样，在训练过程中，损失并不是平滑的。

接下来，用训练好的模型进行预测：

x_data = x.data.numpy() # 获得x包裹的数据
plt.figure(figsize = (10, 7)) # 设定绘图窗口大小
xplot, = plt.plot(x_data, y.data.numpy(), 'o') # 绘制原始数据
yplot, = plt.plot(x_data, predictions.data.numpy()) # 绘制拟合数据
plt.xlabel('X') # 更改坐标轴标注
plt.ylabel('Y') # 更改坐标轴标注
plt.legend([xplot, yplot],['Data', 'Prediction under 1000000 epochs']) # 绘制图例
plt.show()

可以看到，在开头的训练集部分表现良好，但后面就不行了。

书上说，这是没有做归一化导致的：

我们知道，x的取值范围是1~50，而所有权重和偏置的初始值都是设定在(-1, 1)的正态分布随机数，那么输入层到隐含层节点的数值范围就成了的-50~50，要想将sigmoid函数的多个峰值调节到我们期望的位置，需要耗费很多计算时间。事实上，如果让训练时间更长些，我们可以将曲线后面的部分拟合得很好。

重新修改输入层，做归一化：

x = torch.FloatTensor(np.arange(len(counts), dtype = float) / len(counts))

重新进行训练，损失曲线：

loss: tensor(2238.6128, grad_fn=<MeanBackward0>)
loss: tensor(942.9857, grad_fn=<MeanBackward0>)
loss: tensor(713.7626, grad_fn=<MeanBackward0>)
loss: tensor(524.3998, grad_fn=<MeanBackward0>)
loss: tensor(439.8004, grad_fn=<MeanBackward0>)
loss: tensor(235.2034, grad_fn=<MeanBackward0>)
loss: tensor(129.9042, grad_fn=<MeanBackward0>)
loss: tensor(80.1204, grad_fn=<MeanBackward0>)
loss: tensor(56.3744, grad_fn=<MeanBackward0>)
loss: tensor(47.6922, grad_fn=<MeanBackward0>)

在相同训练迭代下，实现了更好的收敛效果。再次绘制预测曲线，这次能够比较好的拟合了：

但这只是在训练集上做“预测”。接下来是实际的预测，即在验证集上预测。采取的验证方法，是从数据集中取出随后 50 条数据，做同样的归一化处理：

# 读取待预测的后面50个数据点
counts_predict = rides['cnt'][50:100] 
x = torch.FloatTensor((np.arange(len(counts_predict), dtype = float) + len(counts)) / len(counts))

# 读取后面50个点的y数值，不需要做归一化
y = torch.FloatTensor(np.array(counts_predict, dtype = float))

# 用x预测y
hidden = x.expand(sz, len(x)).t() * weights.expand(len(x), sz) # 从输入层到隐含层的计算
hidden = torch.sigmoid(hidden) # 将sigmoid函数应用在隐含层的每一个神经元上
predictions = hidden.mm(weights2) # 从隐含层输出到输出层，计算得到最终预测值
loss = torch.mean((predictions - y) ** 2) # 计算预测数据上的损失函数
print(loss)

# 将预测曲线绘制出来
x_data = x.data.numpy() # 获得x包裹的数据
plt.figure(figsize = (10, 7)) # 设定绘图窗口大小
xplot, = plt.plot(x_data, y.data.numpy(), 'o') # 绘制原始数据
yplot, = plt.plot(x_data, predictions.data.numpy()) # 绘制拟合数据
plt.xlabel('X') # 更改坐标轴标注
plt.ylabel('Y') # 更改坐标轴标注
plt.legend([xplot, yplot],['Data', 'Prediction']) # 绘制图例
plt.show()

可以看到，预测完全失败！

过拟合

原因是过拟合了。书上说：

所谓过拟合（overfitting）现象，是指模型可以在训练数据上进行非常好的预测，但在全新的测试数据上表现不佳。

在本例中，原因在与特征变量选错了：

原因就在于我们选择了错误的特征变量：我们尝试用数据的下标（1, 2, 3, …）或者它的归一化（0.1, 0.2, …）来对y进行预测。然而曲线的波动模式（也就是单车的使用数量）显然并不依赖于下标，而是依赖于诸如天气、风速、星期几和是否是节假日等因素。然而，我们不管三七二十一，硬要用强大的人工神经网络来拟合整条曲线，这自然就导致了过拟合的现象，而且是非常严重的过拟合。

单车预测器 2.0

用数据集中的各种环境因素（天气、风速、星期几）来预测。

数据集中的变量分为两种：

类型变量：不连续
数值变量：连续

不能直接用。类型变量转独热编码，将类型转为向量，“独热编码的向量就对应了不同的激活模式。这样的数据更容易被神经网络处理。”。

对数据集进行独热编码转换：

# 所有类型编码变量的名称
dummy_fields = ['season', 'weathersit', 'mnth', 'hr', 'weekday'] 
for each in dummy_fields:
	# 取出所有类型变量，并将它们转变为独热编码
	dummies = pd.get_dummies(rides[each], prefix=each, drop_first=False)
	# 将新的独热编码变量与原有的所有变量合并到一起
	rides = pd.concat([rides, dummies], axis=1)

# 将原来的类型变量从数据表中删除
fields_to_drop = ['instant', 'dteday', 'season', 'weathersit', 'weekday', 'atemp', 'mnth', 'workingday',
 'hr'] # 要删除的类型变量的名称
# 将它们从数据库的变量中删除
data = rides.drop(fields_to_drop, axis=1)

对于数值变量，用变量的均值和标准差进行归一化，针对温度，变成 [-1, 1] 内的分布。

quant_features = ['cnt', 'temp', 'hum', 'windspeed'] # 数值类型变量的名称
scaled_features = {} # 将每一个变量的均值和方差都存储到scaled_features变量中
for each in quant_features:
	# 计算这些变量的均值和方差
	mean, std = data[each].mean(), data[each].std()
	scaled_features[each] = [mean, std]
	# 对每一个变量进行标准化
	data.loc[:, each] = (data[each] - mean)/std

现在的数据集：17397 条数据，59 个变量。最后 504 条数据作为验证集，其它用于训练。

数据集划分：

test_data = data[-21*24:] # 选出训练集
train_data = data[:-21*24] # 选出测试集
# 目标列包含的字段
target_fields = ['cnt','casual', 'registered']
# 将训练集划分成特征变量列和目标特征列
features, targets = train_data.drop(target_fields, axis=1), train_data[target_fields]
# 将测试集划分成特征变量列和目标特征列
test_features, test_targets = test_data.drop(target_fields, axis=1), test_data[target_fields]
# 将数据类型转换为NumPy数组
X = features.values # 将数据从pandas dataframe转换为NumPy
Y = targets['cnt'].values
Y = Y.astype(float)
Y = np.reshape(Y, [len(Y),1])
losses = []

构建神经网络，结构如下：

输入层：56 个神经元
隐含层：10 个神经元
输出层：1 各神经元

使用 PyTorch 自带的神经网络搭建代码：

# 定义神经网络架构，features.shape[1]个输入单元，10个隐含单元，1个输出单元
input_size = features.shape[1]
hidden_size = 10
output_size = 1
batch_size = 128
neu = torch.nn.Sequential(
	torch.nn.Linear(input_size, hidden_size),
	torch.nn.Sigmoid(),
	torch.nn.Linear(hidden_size, output_size),
)

使用 PyTorch 自带的损失函数：

cost = torch.nn.MSELoss()

使用 PyTorch 自带的 SGD 优化器：

optimizer = torch.optim.SGD(neu.parameters(), lr = 0.01)

数据批处理

但是，现在的数据量是16 875条，在这么大的数据量下，如果在每个训练周期都处理所有数据，则会出现运算速度过慢、迭代可能不收敛等问题。

解决方法通常是采取批处理（batch processing）的模式，也就是将所有的数据记录划分成一个批次大小（Batch Size）的小数据集，然后在每个训练周期给神经网络输入一批数据，如图3.22所示。批次的大小依问题的复杂度和数据量的大小而定，在本例中，我们设定batch_size=128。

训练代码：

# 神经网络训练循环
losses = []
for i in range(1000):
	# 每128个样本点划分为一批，在循环的时候一批一批地读取
	batch_loss = []
	# start和end分别是提取一批数据的起始下标和终止下标
	for start in range(0, len(X), batch_size):
		end = start + batch_size if start + batch_size < len(X) else len(X)
		xx = torch.FloatTensor(X[start:end])
		yy = torch.FloatTensor(Y[start:end])
		predict = neu(xx)
		loss = cost(predict, yy)
		optimizer.zero_grad()
		loss.backward()
		optimizer.step()
		batch_loss.append(loss.data.numpy())
	# 每隔100步输出损失值
	if i % 100==0:
		losses.append(np.mean(batch_loss))
		print(i, np.mean(batch_loss))
# 打印输出损失值
plt.plot(np.arange(len(losses))*100,losses)
plt.xlabel('epoch')
plt.ylabel('MSE')

训练过程：

0 0.924881
100 0.26576263
200 0.24588603
300 0.20510724
400 0.12540665
500 0.08741305
600 0.07506938
700 0.069889754
800 0.06674685
900 0.064443044

误差曲线：

在测试集上进行预测：

targets = test_targets['cnt'] # 读取测试集的cnt数值
targets = targets.values.reshape([len(targets),1]) # 将数据转换成合适的张量形式
targets = targets.astype(float) # 保证数据为实数

x = torch.FloatTensor(test_features.values)
y = torch.FloatTensor(targets)

# 用神经网络进行预测
predict = neu(x)
predict = predict.data.numpy()

fig, ax = plt.subplots(figsize = (10, 7))

mean, std = scaled_features['cnt']
ax.plot(predict * std + mean, label='Prediction')
ax.plot(targets * std + mean, label='Data')
ax.legend()
ax.set_xlabel('Date-time')
ax.set_ylabel('Counts')
dates = pd.to_datetime(rides.loc[test_data.index]['dteday'])
dates = dates.apply(lambda d: d.strftime('%b %d'))
ax.set_xticks24]
_ = ax.set_xticklabels24], rotation=45

预测曲线：

可以看到基本是吻合的。但是后半部分效果变差，根据书上说，这与圣诞节有关。

神经网络剖析

神经网络是作为黑箱进行训练的，从研究角度，需要打开黑箱，研究网络的内部情况。

首先，定义函数用于提取网络中的所有参数：

def feature(X, net):
	# 定义一个函数，用于提取网络的权重信息，
	# 所有的网络参数信息全部存储在Neu的named_parameters集合中
	X = torch.from_numpy(X).type(torch.FloatTensor)
	dic = dict(net.named_parameters()) # 从这个集合中提取数据
	weights = dic['0.weight'] # 可以按照“层数.名称”来索引集合中的相应参数值
	biases = dic['0.bias']
	h = torch.sigmoid(X.mm(weights.t()) + biases.expand([len(X), len(biases)]))
	# 隐含层的计算过程
	return h # 输出层的计算

获取预测中表现不好的 3 天的数据：

bool1 = rides['dteday'] == '2012-12-22'
bool2 = rides['dteday'] == '2012-12-23'
bool3 = rides['dteday'] == '2012-12-24'
# 将3个布尔型数组求与
bools = [any(tup) for tup in zip(bool1,bool2,bool3) ]
# 将相应的变量取出
subset = test_features.loc[rides[bools].index]
subtargets = test_targets.loc[rides[bools].index]
subtargets = subtargets['cnt']
subtargets = subtargets.values.reshape([len(subtargets),1])

将这三天数据输入网络，通过函数提取出隐含神经元中的激活数值，还做了归一化还原，方便阅读：

# 将数据输入到神经网络中，读取隐含神经元的激活数值，存入results中
results = feature(subset.values, neu).data.numpy()
# 这些数据对应的预测值（输出层）
predict = neu(torch.FloatTensor(subset.values)).data.numpy()
# 将预测值还原为原始数据的数值范围
mean, std = scaled_features['cnt']
predict = predict * std + mean
subtargets = subtargets * std + mean

接下来，绘制隐含神经元的激活情况：

# 将所有的神经元激活水平画在同一张图上
fig, ax = plt.subplots(figsize = (8, 6))

# 绘制所有隐含神经元
ax.plot(results[:,:],'.:',alpha = 0.3)

# 绘制预测
ax.plot((predict - min(predict)) / (max(predict) - min(predict)),'bo-',label='Prediction')

# 绘制实际
ax.plot((subtargets - min(predict)) / (max(predict) - min(predict)),'ro-',label='Real')

# 绘制第六个神经元
ax.plot(results[:, 5],'.:',alpha=1,label='Neuro 6')
ax.set_xlim(right=len(predict))
ax.legend()
plt.ylabel('Normalized Values')
dates = pd.to_datetime(rides.loc[subset.index]['dteday'])
dates = dates.apply(lambda d: d.strftime('%b %d'))
ax.set_xticks24]
_ = ax.set_xticklabels24], rotation=45

Neuro 6

在书中的训练结果中，Neuro 6 与实际输出曲线比较接近。

但是在我的训练版本中，Neuro 6 始终接近 1，该神经元处于高度兴奋状态。

为什么跟书中不一样呢？

这说明，即使是同样的训练数据、同样的神经网络结构、同样的训练结果，但不同初始状态选练出的模型内部参数也是不同的。

除非我采用跟书中相同的随机数种子、同样的初始状态，才能完全复现书中结果。

下面总结书中分析 Neuro 6 的思路：

对 Neuro 6，绘制不同特征所对应的权重q
横轴代表不同的权重，也就是输入神经元的编号；纵轴代表神经网络训练后的连边权重。
横轴表示各个输入特征
纵轴高或者低表示神经元被激活，表示特征生效

根据分析，作者推断：

与此相对的是，神经元在weathersit_3和hr_6这两个输入上的权重值为负值，并且刚好是低谷，这意味着该神经元会在下雨或下雪，以及早上6点的时候被抑制。通过翻看日历可知，2012年的12月22日和23日刚好是周六和周日，因此Neuro 6被激活了，它们对正确预测这两天的正午高峰做了贡献。

但是，由于圣诞节即将到来，人们可能早早回去为圣诞做准备，因此这个周末比较特殊，并未出现往常周末的大量骑行需求，于是Neuro 6给出的激活值导致了正午单车预测值过高。

通过对单个神经元的分析，能够看出它关注输入属性的哪些部分。

本文作者：Maeiee

本文链接：PyTorch 单车预测器

版权声明：如无特别声明，本文即为原创文章，版权归 Maeiee 所有，未经允许不得转载！

喜欢我文章的朋友请随缘打赏，鼓励我创作更多更好的作品！